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Abstract: 

We  introduce  a  new  signal  model,  called  (K,  G')-sparse,  to 
capture  AT -sparse  signals  in  N  dimensions  whose  nonzero 
coefficients  are  contained  within  at  most  C  clusters,  with 
C  <  K  <C  N.  In  contrast  to  the  existing  work  in 
the  sparse  approximation  and  compressive  sensing  liter¬ 
ature  on  block  sparsity,  no  prior  knowledge  of  the  loca¬ 
tions  and  sizes  of  the  clusters  is  assumed.  We  prove  that 
O  (K  +  CTog(.ZV/C)))  random  projections  are  sufficient 
for  (AT,  C)  -model  sparse  signal  recovery  based  on  sub¬ 
space  enumeration.  We  also  provide  a  robust  polynomial¬ 
time  recovery  algorithm  for  (AT,  C)-model  sparse  signals 
with  provable  estimation  guarantees. 

1.  Introduction 

Compressive  sensing  (CS)  is  an  alternative  to  Shan- 
non/Nyquist  sampling  for  the  acquisition  of  sparse  or 
compressible  signals  in  an  appropriate  basis  [1,2],  By 
sparse,  we  mean  that  only  AT  of  the  N  basis  coefficients 
are  nonzero,  where  K  <C  N.  By  compressible,  we  mean 
the  basis  coefficients,  when  sorted,  decay  rapidly  enough 
to  zero  so  that  they  can  be  well-approximated  as  AT-sparse. 
Instead  of  taking  periodic  samples  of  a  signal,  CS  mea¬ 
sures  inner  products  with  random  vectors  and  then  recov¬ 
ers  the  signal  via  a  sparsity-seeking  convex  optimization 
or  greedy  algorithm.  The  number  of  compressive  mea¬ 
surements  M  necessary  to  recover  a  sparse  signal  under 
this  framework  grows  as  M  =  O  (K  \og{N / K)) 

In  many  applications,  including  imaging  systems  and 
high-speed  analog-to-digital  converters,  such  a  saving  can 
be  dramatic;  however,  the  dimensionality  reduction  from 
N  to  M  is  still  not  on  par  with  state-of-the-art  transform 
coding  systems.  While  many  natural  and  manmade  signals 
can  be  described  to  a  first-order  as  sparse  or  compress¬ 
ible,  their  sparse  supports  often  have  an  underlying  do¬ 
main  specific  structure  [3-6],  Exploiting  this  structure  in 
CS  recovery  has  two  immediate  benefits.  First,  the  number 
of  compressive  measurements  required  for  stable  recovery 
decreases  due  to  the  reduction  in  the  degrees  of  freedom  of 
a  sparse  or  compressible  signal.  Second,  true  signal  infor¬ 
mation  can  be  better  differentiated  from  recovery  artifacts 
during  signal  recovery,  which  increases  recovery  robust¬ 
ness.  Only  by  exploiting  a  priori  information  on  coeffi¬ 
cient  structure  in  addition  to  signal  sparsity,  can  CS  hope 
to  be  competitive  with  the  state-of-the-art  transform  cod¬ 


ing  algorithms  for  dimensionality  reduction. 

Fortunately,  it  is  possible  to  design  CS  recovery  al¬ 
gorithms  that  exploit  the  knowledge  of  structured  spar¬ 
sity  models  with  provable  performance  guarantees  [3, 5, 
6],  In  particular,  the  model-based  CS  recovery  framework 
in  [3]  generalizes  to  any  structured-sparsity  model  that  has 
a  tractable  model-based  approximation  algorithm.  This 
framework  has  been  applied  productively  to  two  struc¬ 
tured  signal  models:  block  sparsity  and  wavelet  trees 
with  robust  recovery  guarantees  from  O  (AT)  measure¬ 
ments  [3],  To  recover  signals  that  have  structured  spar¬ 
sity,  problem-specific  convex  relaxation  approaches  are 
also  used  in  the  literature  with  recovery  guarantees  sim¬ 
ilar  to  those  in  [3];  e.g.,  for  block  sparse  signals,  see  [5, 6], 

In  this  paper,  we  introduce  a  new  structured  sparsity 
model,  called  the  (AT,  C)-model,  that  constrains  the  AT- 
sparse  signal  coefficients  to  be  contained  within  at  most 
C'-clusters.  In  contrast  to  the  block  sparsity  model  in  [5, 
6],  our  proposed  model  does  not  assume  prior  knowledge 
of  the  locations  and  sizes  of  the  coefficient  clusters.  We 
show  that  O  ( K  +  C\og(N/C)))  random  projections  are 
sufficient  for  ( K ,  C)-model  signal  recovery  using  a  sub¬ 
space  counting  argument.  We  also  provide  a  polynomial¬ 
time  model-based  approximation  algorithm  based  on  dy¬ 
namic  programming  and  a  CS  recovery  algorithm  based 
on  the  model-based  recovery  framework  of  [3],  In  con¬ 
trast  to  the  clustered  sparse  recovery  algorithm  based  on 
the  probabilistic  Ising  model  in  [7],  the  (AT,  C)-model  has 
provable  performance  guarantees. 

The  paper  is  organized  as  follows.  Section  2  provides 
the  necessary  theoretical  and  algorithmic  background  on 
model-based  CS.  Section  3  introduces  the  (AT,  C)-model, 
derives  its  sampling  bound  for  CS  recovery,  and  describes 
a  dynamic  programming  solution  for  optimal  (A',  C)- 
model  approximation.  Section  4  discusses  the  aspect  of 
compressibility  and  highlights  some  connections  to  the 
block  sparsity  model.  Simulation  results  are  given  in 
Section  5  to  demonstrate  the  effectiveness  of  the  (A',  C)- 
model.  Section  6  provides  our  conclusions. 

2.  Model-based  CS  Background 

A  AT-sparse  signal  vector  x  lives  in  E k  C  R^,  which 
is  a  union  of  (^)  subspaces  of  dimension  AT.  Other  than 
its  A’-sparsity,  there  are  no  further  constraints  on  the  sup¬ 
port  or  values  of  its  coefficients.  A  union-of-subspaces 


signal  model  (a  signal  model  in  the  sequel  for  brevity)  en¬ 
dows  the  I\  -sparse  signal  x  with  additional  structure  that 
allows  certain  A'-dimensional  subspaces  in  Yjk  and  disal¬ 
lows  others  [4, 8], 

More  formally,  let  x\n  represent  the  entries  of  x  cor¬ 
responding  to  the  set  of  indices  O  C  {1, . . . ,  N },  and  let 

denote  the  complement  of  the  set  fi.  A  signal  model 
A4k  is  then  defined  as  the  union  of  uik  canonical  K- 
dimensional  subspaces 
mK 

Mk  =  |J  Xm,  Xm  :=  {x  :  x\nm  G  Rk,;e|qc  =  0}. 

m—1 

Each  subspace  Xm  contains  all  signals  x  with  supp(a:)  G 
fim.  Thus,  the  signal  model  A4k  is  defined  by  the  set  of 
possible  supports  {Qi, . . . ,  flmK }.  Signals  from  A4k  are 
called  K-model  sparse.  Likewise,  we  may  define  Xi'K 
to  be  the  set  of  c-wise  differences  of  signals  belonging 
to  Mk-  Clearly,  M.k  C  and  AfJj-  C  E^.  In 
the  sequel,  we  will  use  an  algorithm  M(a;;  AT )  that  returns 
the  best  AT-term  approximation  of  the  signal  x  under  the 
model  Mk- 

If  we  know  that  the  signal  x  being  acquired  is  K- 
model  sparse,  then  we  can  relax  the  standard  restricted 
isometry  property  (RIP)  [1]  of  the  CS  measurement  matrix 
$  and  still  achieve  stable  recovery  from  the  compressive 
measurements  y  =  <1>:/;.  The  model-based  RIP  A4  a- -RIP 
requires  that 

(1  -  6mk) IMli  <  \\$x\\l  <  (1  +  SMk)\\x\\1  (1) 

hold  for  signals  x  G  A4 k  [4, 8],  where  5mk  Is  the  model- 
based  RIP  constant. 

Blumensath  and  Davies  [4]  have  quantified  the  num¬ 
ber  of  measurements  M  necessary  for  a  subgaussian  CS 
matrix  to  have  the  Adje-RIP  with  constant  8mk  and  with 
probability  1  —  e-t  to  be 

M  >  —p; —  ( ln(2  ttik)  +  K  In  - - b  t]  .  (2) 

cSMk  V  <W  / 

This  bound  can  be  used  to  recover  the  conventional  CS 
result  by  substituting  mx  =  (^)  ~  ( Ne/K)K . 

To  take  practical  advantage  of  signal  models  in  CS, 
we  can  integrate  them  into  a  standard  CS  recovery  algo¬ 
rithm  based  on  iterative  greedy  approximation.  The  key 
modification  is  surprisingly  simple  [3]:  we  merely  replace 
the  best  AT-term  approximation  step  with  the  best  A-term 
model-based  approximation  M(x;  K).  For  example,  in  the 
CoSaMP  algorithm  [9],  the  best  L  A -term  approximation 
(with  L  a  small  integer)  is  modified  to  incorporate  a  best 
A  AT-term  model -based  approximation.  The  resulting  al¬ 
gorithm  (see  [3])  then  inherits  the  following  model-based 
CS  recovery  guarantee  at  each  iteration  i,  when  the  mea¬ 
surement  matrix  $  has  the  Af^-RIP  with  SM4_  <  0.1: 

||z  -  Xih  <  2— * || at|| 2  +  20^||x  -  xmkWi 

+  -^\\x  ~  XmAi  +  IKIb 


where  xmk  —  M(a:;  K )  is  the  best  model-based  approxi¬ 
mation  of  x  within  A in- 


3.  The  (AT,  C) -Model 


Motivation:  The  block  sparsity  model  is  used  in  appli¬ 
cations  where  the  significant  coefficients  of  a  sparse  signal 
appear  in  designated  blocks  on  the  ambient  signal  dimen¬ 
sion,  e.g.,  group  sparse  regression  problems,  DNA  mi¬ 
croarrays,  MIMO  channel  equalization,  source  localiza¬ 
tion  in  sensor  networks,  and  magnetoencephalography  [3, 
5,6, 10-14],  It  has  been  shown  that  recovery  algorithms 
provably  improve  standard  CS  recovery  by  exploiting  this 
block-sparse  structure  [3,5]. 

The  (AT,  C')-model  generalizes  the  block  sparsity 
model  by  allowing  the  significant  coefficients  of  a  sparse 
signal  to  appear  in  at  most  C  clusters  of  unknown  size  and 
location  (Figure  1(a)).  This  way,  the  ( K .  G'j-modcl  fur¬ 
ther  accommodates  additional  applications  in,  e.g.,  neuro¬ 
science  problems  that  are  involved  with  decoding  of  nat¬ 
ural  images  in  the  primary  visual  cortex  (VI)  or  under¬ 
standing  the  statistical  behavior  of  groups  of  neurons  in 
the  retina  [15].  In  this  section,  we  formulate  the  (AT,  C )- 
model  as  a  union  of  subspaces  and  pose  an  approximation 
algorithm  on  this  union  of  subspaces. 

To  define  the  set  of  (AT,  C)-sparse  signals,  without 
loss  of  generality,  we  focus  on  canonically  sparse  signals 
in  AT  +  2  dimensions  whose  first  and  last  coefficients  are 
zero.  Consider  expressing  the  support  of  such  signals  via 
run-length  coding  with  a  vector  (3  =  (pi , . . .  ,/?2C+i) 
(Pj  0),  where  /30dd  counts  the  number  of  continuous 
zero-signal  values  and  peven  counts  the  number  of  contin¬ 
uous  nonzero-signal  values  (i.e.,  clusters). 

Definition:  The  ( K,  C)-sparse  signal  model  M(k,c)  Is 

defined  as 


M 


( K,C ) 


=  {  x  £  R 


N+2 


2C+1 


Y.  Pi  =  N  +  2,Yfoi  =  K  >■ 


(3) 


Sampling  Bound:  The  number  of  subspaces  rri(K,c)  in 
M  (k,c)  can  t*e  obtained  by  counting  the  number  of  posi¬ 
tive  solutions  to  the  following  integer  equations: 


Pi+  @2  +  ■•■  +  P2C+1  —  N  +  2, 
P2  +  Pi  +  •  •  ■  +  P2C  =  K, 


which  can  be  rewritten  as 

01+ fo  + ■■■  +  P2C+1  =  N  +  2  —  K, 

P2  +  Pi  +  •  •  •  T  P2C  =  A . 

Note  that  the  number  of  positive  integer  solutions  to  the 
following  problem: 

Pi  +  02  +  03  +  •  ■  •  +  Pn  =  AT, 

is  given  by  (^Ci)  ■  Then,  we  can  count  the  solutions  to  the 
two  of  decoupled  problems  in  (4)  and  multiply  the  number 
of  solutions  to  obtain  rri(K,c): 


/AT  +  1  -  K\  /AT  -  1\ 

V  ){c-ij- 


m(K,C)  — 


(5) 


Plugging  (5)  into  (2),  we  obtain  the  sampling  bound 
for  M{k,c)- 


M  =  O 


K  +  C  log 


N 

C 


(6) 


Note  that  the  ( K ,  C)-sampling  bound  (6)  becomes  the 
standard  CS  bound  of  M  =  O  (if  log  when  C  ~  K . 
Model  Approximation  Algorithm:  In  this  section  we 
focus  on  designing  an  algorithm  Ml  fa:;  K ,  C)  for  finding 
the  best  (if,  C)-model  approximation  to  a  given  signal 
x.  The  algorithm  uses  the  principle  of  dynamic  program¬ 
ming  [16].  For  simplicity,  we  focus  on  the  problem  of 
finding  the  cost  of  the  best  (if,  C)-clustered  signal  ap¬ 
proximation  in  l2.  This  solution  generalizes  to  the  best 
(if,  C')-clustered  signal  approximation  in  £p  for  p  >  1. 
The  actual  sparsity  pattern  can  be  then  recovered  using 
standard  back-tracing  techniques;  see  [16]  for  the  details. 

The  algorithm  M(x;  K ,  C)  computes  an  array 
cost  [i,j,k,c],  where  l<i<j<N,0<k<K, 
and  0  <  c  <  C.  At  the  end  of  the  algorithm,  each  entry 
cost[i,  j,  k,  c]  contains  the  smallest  cost  of  approximating 
Xi-j,  the  signal  vector  restricted  to  the  index  set  [i,. . . ,  j], 
using  at  most  k  non-zero  entries  that  span  at  most  c  clus¬ 
ters.  M(x;  if,  C)  performs  the  following  operations. 

( Initialization )  When  either  c  =  0  or  k  =  0,  the  sig¬ 
nal  approximation  costs  can  be  computed  directly,  since 
cost[i,j,0,c]  =  \\xi:j\\l  and  cost[i,  j,  k,  0]  =  ||xj:j|||, 
for  all  valid  indices  i,j,k,c.  Moreover,  for  all  entries 
i.j,  k,  c  such  that  c  >  0  and  j  —  i  +  1  <  fc,  we  have 
cost[i,  j,  k,c]  =0  since  we  can  include  all  j  —  i  +  1  coor¬ 
dinates  of  the  vector  xl:j  in  the  approximation. 


4.  Additional  Remarks 

Compressibility:  Just  as  compressible  signals  are 
nearly  K -sparse  and  live  close  to  the  union  of  sub¬ 
spaces  Y>k  in  M.N ,  (if,  C')-compressible  signals  are  nearly 
( K ,  C)-model  sparse  and  live  close  to  the  restricted  union 
of  subspaces  Ai  (k,c)-  Here,  we  rigorously  introduce  a 
( K ,  C)-compressible  signal  model  in  terms  of  the  decay 
of  their  (if,  C')-model  approximation  error. 

We  first  define  the  I2  error  incurred  by  approximating 
x  £  by  the  best  approximation  in  Ai (k,c): 

<?M(k  C )  ( x )  -  inf  \\x~xh  =  ||z  — M(a:;  if,  C)\\2. 

x£M(k,c ) 

The  decay  of  the  (K,  C)-model  approximation  error  in  (7) 
defines  the  (if,  (^-compressibility  of  a  signal.  Then,  a  set 
of  ( K ,  C)-model  s-compressible  signals  is  given  by 


ms  =  l  x  e  RN  :  <rMj(x,o)  (*)  <  S(jif)“1/S, 

^  s  (7) 

1  <  if  <  N,  S  <  00,  j  =  1,...,  ^7  i. 

Define  Sm  as  the  smallest  value  of  S  for  which  this  con¬ 
dition  holds  for  x  and  s. 

We  use  the  restricted  amplification  property  (RAmP) 
and  the  nested  approximation  property  (NAP)  in  [3]  to 
ensure  that  the  (if,  C)-model  based  CoSaMP  recovery 
possesses  the  following  guarantee  for  ( K ,  C)-model  s- 
compressible  signals  at  each  iteration  i: 


( Main  loop )  All  other  cost  entries  can  then  be  com¬ 
puted  using  the  following  recursion: 

costfi,  j,  k,  d  =  min  min  min 

c*— 0...C  fc*=0...fc  j*=i...j- 1 


cost[i,  j*,  k*,  c*]  x  cost[j*  +  1  ,j,  k  —  k*,  c  —  c*]  > . 


II*  -  Xi || 2  <  2  ‘11*112  +  35  ^||n||2  +  ^-(1  +  In \N/ if]) j  , 

(8) 

when  $  has  the  Ai  ^n-RIP  with  5 <0.1  and  the 
(err,  r)-RAmP  with  ck  <  0.1  and  r  =  s  —  1. 


The  correctness  of  the  algorithm  follows  from  the  follow¬ 
ing  observation.  Let  v  be  the  best  ( k ,  c)-clustered  approx¬ 
imation  of  Xi-.j.  Unless  all  entries  of  xl:j  can  be  included 
in  the  approximation  v  (in  which  case  j  —  i  +  1  >  k  and 
the  entry  has  been  already  computed  during  initialization), 
then  there  must  exist  an  index  l  €  such  that  Xi  is 

not  included  in  v.  Let  l*  =  l  if  l  <  j,  and  l*  =  j  —  1  oth¬ 
erwise.  Let  k*  be  the  number  of  non-zero  entries  present 
in  the  left  segment  of  Vu*,  and  let  c*  be  the  number  of 
clusters  present  in  that  left  segment.  Then,  it  must  be  the 
case  that  Vi-j*  is  the  best  (fc*,  c*) -approximation  to  xp.i, 
and  vi+i-j  is  the  best  (k  —  k*,  c  —  c*)-approximation  to 
X(i*+i):j-  Otherwise,  those  better  approximations  could 
have  been  concatenated  together  to  yield  an  even  better 
(k,  c) -approximation  of  xt:j .  Thus,  the  recursive  formula 
will  identify  the  optimal  split  and  compute  the  optimal  ap¬ 
proximation  cost. 

The  cost  table  contains  O  (N1  KC)  entries.  Each  en¬ 
try  can  be  computed  in  O  (NKC)  time.  Thus,  the  running 
time  of  the  algorithm  is  O  (iV3A'2C2). 


Simulation  via  Block  Sparsity:  It  is  possible  to  recover 
(if,  C)-sparse  signals  by  using  the  block  sparsity  model  if 
we  are  willing  to  pay  an  added  penalty  in  terms  of  the 
number  of  measurements.  To  demonstrate  this,  we  define 
uniform  blocks  of  size  K/C  (e.g.,  average  cluster  length) 
on  the  signal  space.  Then,  it  is  straightforward  to  see  that 
the  number  of  active  blocks  B  in  the  block  sparse  model 
is  upper-bounded  by 

B  <  2(C  —  1)  +  A  X)  <  3 C.  (9) 

To  reach  this  upper  bound,  we  first  construct  a  ( K ,  C)- 
sparse  signal  that  has  (C  —  l)-clusters  with  2  coefficients 
and  a  single  cluster  with  the  remaining  sparse  coefficients. 
We  then  place  the  clusters  with  two  coefficients  at  the 
boundary  of  the  block  sparse  model  so  that  each  cluster  ac¬ 
tivate  two  blocks  in  the  block  sparse  model  to  arrive  at  (9). 
Then,  the  ( K ,  C) -equivalent  block  sparse  model  requires 
M  =  O  ( BK/C  +  B  log  jj)  samples,  where  B  =  3 C. 


Figure  1 :  Monte  Carlo  simulation  results  for  ( K ,  C)-model  based  recovery  with  K  =  10,  C  =  2.. 


5.  Experiments 

In  this  section  we  demonstrate  the  performance  of  (K.  C)- 
model  based  recovery.  Our  test  signals  are  the  class  of 
length-100  clustered-sparse  signals  with  I\  =  10,  C  =  2. 
We  run  both  the  CoSaMP  algorithm  as  well  as  (K .  C)- 
model  based  CoSaMP  algorithm  [3]  until  convergence  for 
1000  independent  trials.  In  Fig.  1(a),  a  sample  realization 
of  the  signal  is  displayed.  It  is  evident  from  Figs.  1(b)  and 
(c)  that  enforcing  the  structured  sparsity  model  in  the  re¬ 
covery  process  significantly  improves  CS  reconstruction 
performance.  In  particular.  Fig.  1(b)  demonstrates  that  ap¬ 
proximately  85%  of  the  signals  are  almost  perfectly  re¬ 
covered  at  M  =  2.5 I\,  whereas  CoSaMP  fails  to  recover 
any  signals  at  this  level  of  measurements.  Instead,  tra¬ 
ditional  sparsity-based  recovery  requires  M  >  4.5 K  to 
attain  comparable  performance.  Similarly,  Figure.  1(c) 
displays  the  rapid  decrease  in  average  recovery  distortion 
of  our  proposed  method,  as  compared  to  the  conventional 
approach.  The  {K ,  C)-sparse  approximation  algorithm 
codes  are  available  at  dsp. rice.edu/software/KC. 

6.  Conclusions 

In  this  paper,  we  have  introduced  a  new  sparse  signal 
model  that  generalizes  the  block-sparsity  model  used  in 
the  CS  literature.  To  exploit  the  provable  model-based  CS 
recovery  framework  of  [3],  we  developed  a  dynamic  pro¬ 
gramming  algorithm  that  computes,  for  any  given  signal, 
its  optimal  /^-approximation  within  our  clustered  spar¬ 
sity  model.  We  then  demonstrated  that  significant  per¬ 
formance  gains  can  be  made  by  exploiting  the  clustered 
signal  model  beyond  the  simplistic  sparse  model  that  are 
prevalent  the  CS  literature. 
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